Demixing behavior in two-dimensional mixtures of anisotropic hard bodies 
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Scaled particle theory for a binary mixture of hard discorect angles and for a binary mixture 
of hard rectangles is used to predict possible liquid-crystal demixing scenarios in two dimensions. 
Through a bifurcation analysis from the isotropic phase, it is shown that isotropic-nematic demixing 
is possible in two-dimensional liquid-crystal mixtures composed of hard convex bodies. This bifur- 
cation analysis is tested against exact calculations of the phase diagrams in the framework of the 
restricted-orientation two-dimensional model (Zwanzig model). Phase diagrams of a binary mixture 
of hard discorectangles are calculated through the parametrization of the orientational distribution 
functions. The results show not only isotropic-nematic, but also nematic-nematic demixing ending 
in a critical point, as well as an isotropic-nematic-nematic triple point for a mixture of hard disks 
and hard discorectangles. 



PACS numbers: 64.70.Md,64.75.+g,61.20.Gy 



I. INTRODUCTION 

The demixing behavior of hard-core three-dimensional 
additive mixtures composed of particles without orienta- 
tional degrees of freedom, in particular the hard-sphere 
and parallel hard-cube systems, has been recently stud- 
ied in depth by several authors using Monte Carlo sim- 
ulation techniques Q, 0, and theory 0, |j| . The main 
conclusions that can be drawn from these studies is that 
fluid-fluid demixing is always mctastablc with respect to 
fluid-solid demixing for asymmetric mixtures in which 
the solid phase is composed of big particles and the fluid 
phase is enriched in small particles (demixing behavior 
similar to the one found in mixtures of parallel hard cubes 
on a cubic lattice [(|). The physics behind this demix- 
ing behavior is known as the depletion effect 0, and 
can be explained as the effective attraction between two 
large particles due to the uncompensated osmotic pres- 
sure exerted by the small particles when the two excluded 
volumes between big and small particles overlap. Thus 
in the demixed phases enriched by each of the species 
the accessible volumes to the particles is maximized and 
as a consequence the total configurational entropy is in- 
creased. 

Three-dimensional mixtures of additive anisotropic 
particles possess a demixing behavior which differs from 
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that of mixtures composed of isotropic particles because 
of the presence of particle orientational degrees of free- 
dom, ft has been shown that demixing strongly depends 
on the shape of the particles (spherical, oblate or prolate) 
II H Elll Elll 13 and for a particular geometry, 
on its aspect ratio (the ratio between the characteris- 
tic lengths of particles) 0, 0, 0|. For example, one 
of the demixed phases can be an oriented phase, i.e., a 
nematic (N) phase where the particles are aligned on av- 
erage along the nematic director, ft is easy to predict 
the most common demixing scenario, it will proceed be- 
tween a N phase composed of the longer particles and an 
isotropic (I) phase composed of the short-particle compo- 
nent, fn this case the excluded volume between the large 
particles is minimized by their parallel alignment. I-I and 
N-N demixing was also found in mixtures of anisotropic 
particles 0IIMII1G2II1I3II1IIIII3, 

showing that 

fluid-fluid demixing is a common scenario in these mix- 
tures, which has recently been confirmed experimentally 



The scaled particle theory (SPT) was initially de- 
veloped for hard spheres If9l and later extended to 
anisotropic particles [2TI l'2ll l22l | . The usual formulation 
of SPT for mixtures of hard convex bodies restricted to 
the isotropic orientational phase has as a main ingredi- 
ent the expression for the second virial coefficients as a 
function of the volumes, surface areas, and mean curva- 
tures of hard convex bodies [2^, 0] ■ This exact result 
obtained in the 1950's [23| was used to show that in two 
dimensions f-I demixing is not allowed |25| . the gain in 
accessible volume after demixing is much lower in two 



2 



dimensions, which explains the stability of mixtures with 
respect to I-I demixing at any composition. This result 
has been confirmed by other theories constructed from 
the expression of the direct correlation function in terms 
of the geometric measures of the particles [lOj, in the 
same spirit as the fundamental measure theory (FMT) 
for hard spheres [26l |27| . Some authors have studied the 
possibility of demixing in a mixture of p erfectly oriented 
two-dimensional additive particles [lOj, arriving at the 
same conclusion, the mixture is always stable. However, 
in contrast to these findings, we will show in section III 
that the SPT approximation does predict N-N demix- 
ing in the limit of parallel alignment. Also, the previous 
negative results cannot be taken as definite, since the 
suppression of orientational fluctuations is a crude ap- 
proximation to study the possibility of N-N demixing, so 
that the question about the existence of demixing (I-N 
or N-N) in two dimensions is still open. Note that we 
explicitly distinguish I-N demixing from the usual orien- 
tational I-N phase transition. The reason for this will be 
explained later in Sec. IV C. 

As opposed to additive mixtures, nonadditive mixtures 
of hard disks can demix, as was shown by several authors 
[28l |2^| . This result is not surprising since the demixed 
phases, composed of practically the same kind of species, 
are approximately additive and, as a consequence, have 
less excluded volumes between any pair of particles com- 
pared to the mixed state. Nonadditive mixtures made 
up of anisotropic particles have not been analyzed so far. 
The coupling between nonadditiveness in the interactions 
and orientational degrees of freedom is expected to give 
rise to interesting phenomenology, which is worth inves- 
tigating. We will not pursue these aspects in the present 
paper. 

The purpose of this work is to shed light on the ques- 
tion about the existence of I-N or N-N demixing in two- 
dimensional mixtures of additive hard anisotropic par- 
ticles. We have used the SPT formalism specified for 
a mixture of hard rectangles (HR) and hard discorect- 
angles (HDR), and applied a bifurcation analysis from 
the isotropic phase to find the possibility of I-N demix- 
ing. This analysis has been shown to be useful in the 
study of phase transitions with symmetry breaking. In 
particular it was used to study the I-N transition in two 
dimensions in a one-component system of hard needles 
|3(J, and in a system composed of self-assembled rods 
|31| in the same two-dimensional Onsager limit. Using 
the bifurcation analysis, we have found that, even for a 
system of HDR (where the I-N transition is continuous 
in the one-component limit, as indicated by simulation 
|32j and density- functional theory |33|) the mixture can 
demix in two phases of different composition. In addi- 
tion, we have calculated the phase diagrams of different 
mixtures of HDR and confirmed this demixing scenario. 
In some mixtures we have also found a N-N phase sepa- 
ration which ends in a critical point, as well as a triple 
coexistence between an I phase and two different nemat- 



The paper is organized as follows. In Sec. II we present 
the theoretical model, and in Sec. Ill the bifurcation 
analysis. The results from this analysis are shown for 
HDR in Sec. IV A, for HR in Sec. IV B, and a check 
using a very simple model (Zwanzig approximation) is 
presented in Sec. IV C. Section V shows the results from 
the calculations of the phase diagrams of HDR mixtures. 
Finally some conclusions are drawn in Sec. VI. 



II. MODEL 

The key quantity in the development of the SPT for 
a mixture of hard convex bodies in two dimensions (see 
Ref.(34|) is the averaged (over all possible orientations of 
a fluid particle of species v) excluded area between an 
inserted scaled particle s with orientation (f>i (measured 
from the nematic director) and a fluid particle of species 
v with orientation 02, i.e., 

(KV.H^^i) = / d^KihK^JLs,^,^), (l) 



where L s and a s are the length and width of the scaled 
particle, h v ((j)2) is the orientational distribution function 
of species v and <\>\i = 4>\ — <p2 the relative angle between 
the axes of particles s and v. The superindex a labels the 
nature of the particle, either hard rectangles (a = HR) 
or hard discorectangles (a = HDR). The excluded area 
between two rectangles (s and v) is 

K™. l /( i s I fT s,0i2) = {L U L S + a v a s ) | sin0i 2 | + v v + v s 



(L u <t s + L s a v ) I cos 0i2 1 



(2) 



For hard discorectangles 



Kxcl^^) "*)^) = L U L S \ sin ^12 1 + v v + v s 

7T 

+ L v a s + L s a v + — o v o s , (3) 

where for HR vp = Lpap is the area of species /? 
(/3 = {v,s}) while for HDR v = L fj ap + The 
reversible work required to insert the scaled particle with 
fixed orientation coincides with the excess chemical po- 
tential and, in the limit of small sizes (L s <C L Vl a s -C o v 
for any 1/), it has the following asymptotic form 



/3aw(<M ~ ^ (O) (i s ,cr s ,0i) 



= -In 



1 -^PAV£ci,„)(Ls,<7 a 



,(4) 



where p„ is the density of species v and the sum runs 
over all species. In the opposite limit of large sizes 
(L s 3> L u , a s S> ov) this work coincides with the ther- 
modynamic work required to open a cavity of area v s , 
which is equal to Pv s , where P is the fluid pressure. The 
SPT interpolates between both limits using a Taylor ex- 
pansion of the function fj,^ (L s: er s , <f>i) around the value 
(L s ,(j s ) — (0,0). The second term of this expansion is 
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fixed to Pv s . Finally, all the particle lengths are taken 
to be those of any one of the species, say v, which results 
in 



PfJ>exc,v(<l>l) = - ln(l - T}) 

E T Pr/#2^(^)K ( r 0) (fc) 
1 — T) 



+ /3Pv„, (5) 



where T] = Y] p v v v is the total packing fraction and 

Vur\<t>i2) = V£ cl)V (L T ,a T ,<pi2) -v v - v T . The excess 
chemical potential of species v is the angular average 



/3/W,i/ = J d<j) 1 h^(4>i) [/3fi exCtU ((t>i)] = -ln(l - 77) 



1-T) 



0Pv v , 



(G) 



where ((• • ■)) means the following double angular aver- 
age: 



Integrating the thermodynamic relations 
d/3P ^ dfin cxc.r 

with the use of Eq. © allows us to find 



i £, t ^Pt((k ( t 0) )) 

2 (l-ry) 2 



(3P = 



1 — 77 



(8) 



where p = p v is the total density. Finally, through 
the definition of the pressure (3P = p + J2„ Pv [f3pcxc,u] — 
$ cxc , where $ cxc = (3T cyic /V is the excess part of the free 
energy density in reduced units, the result (JHJ, and Eq. 
© we obtain 



-/3ln(l -7?) + - ■ 9) 

2 1 — 77 



The ideal part of the free energy density of the mixture 



is 



$id =11 ( ln ^ - 1 + J Q d( f> h M) In [irfh,(<l>)]J ,(10) 

where all distribution functions are normalized as 
L d(j)h v (4>) — 1 (note that, in view of the head-tail 
symmetry of the particles, the angle <f> can be restricted 
to the interval [0, 7r]). The functional minimization of 
$ = ^id + $exc with respect to the {h u } allows, as usual, 
to find the equilibrium distribution functions and corre- 
spondingly the equilibrium free energy of the mixture. 

The above theoretical scheme will be used in the fol- 
lowing section to develop a bifurcation analysis of I-N 



demixing for HDR and HR. Also it will be used to cal- 
culate the phase diagram of mixtures. In order to do 
that we need to fix the fluid pressure, which means that 
the composition of one of the species and the total fluid 
density are no longer independent variables. Once the 
independent variable is chosen, the dependent variable 
can be calculated through the constant pressure crite- 
rion. Thus, the adequate thermodynamic potential to 
work with is the Gibbs free energy per particle in re- 
duced units j3g = (<& + j3P) /p. The minimization of $ 
with respect to the order parameters (in the case of the 
nematic phase), the condition of constant pressure and 
the double tangent construction on j3g with respect to 
the composition of one of the species, allows us to calcu- 
late the coexistence condition between different phases. 
Changing the pressure and repeating the above steps we 
have found the phase diagrams of different binary mix- 
tures. 



III. BIFURCATION ANALYSIS 

The usual bifurcation analysis for a one-component 
fluid with symmetry breaking includes (i) an order pa- 
rameter expansion of the free energy around the bifurca- 
tion point; (ii) the calculation of the inverse isothermal 
compressibility of the ordered phase at the same point. 
The order of the phase transition can be elucidated by the 
combined use of both criteria, the sign of the free-energy 
difference between the ordered and disordered phases, al- 
ready minimized with respect to the order parameter, 
and the sign of the isothermal compressibility of the or- 
dered phase, evaluated at bifurcation. If the system ex- 
hibits a tricritical point, its location can be obtained from 
the vanishing of either the first coefficient in the expan- 
sion of the free-energy difference or the isothermal com- 
pressibility, depending on which of them occurs first [35| . 

For binary mixtures fluid-fluid demixing without sym- 
metry breaking is usually calculated from the vanishing of 
the determinant of the matrix with elements d 2 <& / dpidpj 
which means that the stability of the mixture with re- 
spect to volume and composition fluctuations is violated. 
This allows us to obtain the demixing spinodal. How- 
ever, in a symmetry-breaking transition the above matrix 
should be calculated from the minimized free-energy of 
the ordered phase using the order-parameter expansion 
up to the order required. This criterion is equivalent to 
the loss of convexity of the Gibbs free energy per particle 
of the ordered phase with respect to the mixture compo- 
sition at the bifurcation point. Since our aim is the study 
of I-N demixing, which involves an orientational symme- 
try breaking, we will implement the latter scheme. 

A Fourier-series decomposition of the orientational dis- 
tribution functions 



M0) 



1 

7T 



4^ cos(2/c0) 

fe>i 



(11) 



4 



should retain only even harmonics, due to the symmetry 
of the particles studied here [h v ((f>) — h v {-K — </>)]. In the 
neighborhood of the I-N bifurcation point we can assume 
that the Fourier amplitudes {M. } are small. A Taylor 
expansion of the difference in free energy per particle 
Alp = <y3N — <fii (<ys = $/p) between N and I phases, up to 
fourth order, is therefore valid (the order is defined by a 
small bifurcation parameter e, hjf 1 ~ e k ). Up to fourth 
order, the ideal contribution to Aip reads 



id 




whereas the excess contribution is 

A((^)) 



We have defined 



with 



A((^ } )) 
(v) 



1 T$ h (p) h {u) 



(12) 



(13) 



(14) 



(15) 



k>l 



where it is understood that the sum is to be truncated 
at fourth order in the expansion parameter. In the above 
expresions we have also defined 



T (k) 



for HDR, 



[L„ + (-l) fc <r M ] [L v + {-l) k a v } 



(16) 



for HR 



with Xp, — Pfj./p being the molar fraction of species //, 
y = 7y/(l — 77), and (v) — Yl v x v v v the average particle 
area of the mixture. The complete free-energy expansion 
up to fourth order reads 



(4) , (1) 

11 I'n 



j(4) 
'22 



(17) 



The expressions for the affi's are 



7 W 

,(3) 



,(4) 



^V x v T (k) 
(4fc 2 -l)7T 

■I' ft 

32' 



, fe = l,2, (18) 
/i,i/= 1,2, (19) 



where is the Kronecker delta. 

Minimizing Aip with respect to all the amplitudes ex- 
cept one (say hi ), substituting the results for the other 
amplitudes in (|17f) . and neglecting all terms with order 
higher than four, we obtain an effective free energy dif- 
ference as 



A<p = A[h\ 



Blh\ 



.00 



(20) 



with the following explicit expressions for the coefficients 
A and B: 



A 



B 




(21) 



,(2) 



V 



J 2 ) 



(22) 



Note that A and _B depend on x = X2 and y. To pro- 
ceed we must bear in mind that demixing phase transi- 
tions usually imply fractionation (different composition 
x of the coexisting phases) as well as a change in packing 
fraction r\. Therefore, we should expand both the compo- 
sition x and the variable y around the bifurcation point 

x^x*+x^^y, y~y*+y {2) (h^y . (23) 

Inserting the expansions i|23|) in the first derivative of 
(j2Tj|l with respect to hi and equating the result to zero, 
order by order, we obtain the following conditions: 



A* = 0, O hY 



A* x x {2) + AtyW = -25*, O 



h\ 



(i) 



(24) 
(25) 



where /* = f(x*,y*) with (/ = A,B,A x ,A y ), and the 
subindices x, y in Eq. (|25|1 mean the partial derivatives of 
A with respect to x and y, respectively. Solving Eq. 124|) 
gives the packing fraction as a function of composition 
at the bifurcation point, whereas Eq. (|25|l allows to find 
a relation between x' 2 ' and y^ 2 \ respectively. 

Expanding (|20|l up to fourth order around the bifurca- 
tion point and using l|24|) and (|25|l we obtain the energy 
difference as 



Atp = -B* (h\ 



(26) 



5 



Depending on the sign of B* the nematic branch bifur- 
cates below (positive sign) or above (negative sign) the 
isotropic branch. The latter case corresponds to a first- 
order transition. Using 12/il) and 1251) . Eq. 12tj|) can be 
rewritten as 

= <pi - ^ [A* x (x - x*) + A* y (y - y*)]\ (27) 

The stability of the mixture with respect to volume and 
composition fluctuations is guaranteed when 



d 2 <S> 



6> 2 $ 



d 2 $> 
dpidp 2 



> 0. 



(28) 



Fixing the areas of all species to one (in this case the 
demixing criterion will depend on the difference in parti- 
cle shapes and not on their areas) the Eq. Ij28(l , in terms 
of the variables (x,y), can be written as 



d_ ( 2 dcp\ 8hp_ 
dy \ dy J dx 2 



_gvy 

dx dy J 



>0, 



(29) 



which, for the N phase and at the bifurcation point [with 
use of i|27fl] , becomes 



(i + yT 

2(y*fB* 



d_ 

dy 



; dyi 
dy 



(S)V)^ + (S?)V4) 



> o. 



(30) 



It is easy to show that H\{x, rj) > for HDR, freely rotat- 
ing HR and HR in the Zwanzig approximation. This re- 
sult confirms the general wisdom |25fl mentioned above on 
the stability of a two-dimensional mixture of hard bodies 
against I-I fluid demixing, in the framework of the SPT. 
However, fixing the particle orientations to be parallel we 
obtain the condition for Eq. I|28|) to be violated as 



i] > rf = f + x(l — x) (t — t 1 y 



-1/2 



(31) 



with 



K2 — 1 



= < \ Kl-1 



4, for HDR, 



(32) 



for HR, 



and, correspondingly, the position of the critical point in 
the pressure-composition plane is 



f 3(t + t~ 1 )-2 

^c^^, PP C VI = - ; — y, 

2 (t + t- 1 -2) 2 



(33) 



where k„ = (L v + a v ) ja v for HDR while k„ = L v la v 
for HR (remember that these results apply to mixtures 
such that all particle areas have been taken to be unity) . 
This result confirms that there exists N-N demixing in 
the approximation that particles are perfectly aligned. 

The condition — will be used in the next section 
to calculate the I-N demixing region in the plane A — x*, 
where A = K2/K1 and n v is the aspect ratio of species 
v, and it is equivalent to the condition that the Gibbs 
free energy per particle of the N phase have zero curva- 
ture (second derivative with respect to X\ ) at bifurcation, 
thereby determining the I-N demixing tricritical point. 



IV. RESULTS FROM THE BIFURCATION 
ANALYSIS 

This section shows the results from the bifurcation 
analysis, whose formalism was introduced in the pre- 
ceding section. We have divided this section into three 
sections. The first two are devoted to the study of the 
freely-rotating HDR and HR models, respectively, while 
in the third section we discuss the implementation of the 
formalism to the restricted-orientation model (Zwanzig 
model). In this last section results from the calculation 
of the exact phase diagrams of this model are shown. The 
main purpose is to check the results from the bifurcation 
analysis against exact calculations (note that there is no 
need to use any parametrization in the Zwanzig model). 



A. Hard discorectangles 

The solution of Eq. (|2"l}> for HDR allows us to write 
the packing fraction value at the bifurcation point as a 
function of the composition of the mixture, 



1 



= 1 



37T (v) 



(34) 



where we have defined (u n ) = x*u™ for a generic 
quantity u. The expression for B* which, as pointed 
out before, defines the relative positions of the 1 and N 
branches near the bifurcation point, is 

B* = ^i-[z 2 (3 + zi)r 2 -2z 1 z 2 r + zi(3 + z 2 )],(35) 

where r = L 2 /L 1 and z v = x* v L 2 v j{L 2 ) < 1. The 
quadratic polynomial P(r) enclosed by the square brack- 
ets has a discriminant D = — 48z\Z2, which obviously is 
always negative, while P(0) > 0. Thus the coefficient 
B* is always positive and, as a consequence, the nematic 
branch always bifurcates from below with respect to the 
isotropic one. 

We take the areas of all particles to be equal to I, 
which means that the particle length and width, in units 

r 1/2 

of Vfj, , are 



-1/2 



(k^-I)^. (36) 
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FIG. 1: Solutions of — in the x* (composition)-A plane, 
where A = Ks/ki, in logarithmic scale, for HDR mixtures with 
particle areas equal to unity. The curves are shown for dif- 
ferent values of «2- The shaded areas represent the demixing 
regions predicted from < 0. From outside to inside the 
values of K2 are 1, 2, 3, 5, and 00 (the Onsager limit). 
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FIG. 2: Reduced pressure f3Pvi versus composition xi for 
a binary mixture of HDR with m = 5 and k,2 = 2. Filled 
circles indicate the intersections of a vertical line located at 
the corresponding value of A = 2/5 and the solution of = 
in Fig. The continuous line is the portion corresponding 
to the shaded region in the latter figure. 



Inserting H34J) and Ij35(l in (|30|l . the final expression for 
can be written in a particularly compact form, 




where p = 2L + 7rer is the perimeter of the particle. 

The equality to zero of the expression enclosed by 
square brackets in IJS7J, for the specific case L 2 = (or 
K2 = 1) i.e. a binary mixture of hard disks and HDR), 
gives us the analytic solution x\ — lhir / '4(67r+Lf), shown 
in Fig. n]as the curve which encloses all the other curves 
in the A — x* plane (A = K2/K1). Note that the above 
expression is not valid when A = 1 or x\ = (the one- 
component hard disk fluid). In Fig. ^the solutions of 
H£j = for other mixtures, including the Onsager limit 
(i.e. a mixture of hard needles), are shown. It is re- 
markable that, even in this limit, the mixture can demix. 
To calculate these curves the aspect ratio k 2 was fixed, 
while Ki was varied from 1 to 10 4 , and = was 
solved for x*. The curve corresponding to the Onsager 
limit — > 00 (jj, — 1,2) was calculated by first tak- 



ing the corresponding limit of the expression enclosed by 
square brackets in Eq. I|37f) (which depends only on A) 
and solving again for x*. 

We now proceed to explain the physics behind the 
behavior shown in Fig. ^ using the particular values 
A = 2/5, k-2 = 2 as an example. A vertical line located 
at this value of A intersects the curve corresponding to 
«2 = 2 at two points, giving two values of x\ at which 
= 0. In Fig. [3 we plot the reduced pressure as 
a function of x± for these values of A and As can 
be seen from this figure the pressure is a monotonically 
decreasing function of x\, which means that the above- 
mentioned intersection with a higher value of x\ has a 
lower pressure than the other one. This, in turn, implies 
that the first point corresponds to a genuine tricritical 
point where, for the first time, the Gibbs free energy loses 
its convexity with respect to the composition variable. 
The second point also corresponds to a loss of convexity 
of the Gibbs free energy but, as will be discussed later, 
lies inside the two-phase region of a demixing transition. 



B. Hard rectangles 

As shown in Refs. [3i| and |34| the system of HR ex- 
hibits a transition to a phase with fourfold symmetry, 
the so-called tetratic phase (N t ). In this phase the ori- 
entational distribution function has a symmetry under 
rotation by tt/2, h{<j>) = h{n/2 - <f>). The I-N t transi- 
tion was always found to be of second order in the whole 
region of its stability (1 < k < 2.21). Thus the ques- 
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tion naturally arises as to the relative stability of this 
phase in the binary mixture. The above symmetry dic- 
tates that the odd Fourier amplitudes {/i 2 j-i} should be 
equal to zero. Using this constraint and carrying out the 
same expansion for the free-energy difference around the 
bifurcation point, we arrive at the same Eq. H17|). with 
the coefficients affi (k = 1, 2) obtained from Eq. l|18|l . but 
through the substitution k — > 2k on the right-hand side 
of JTSJ, while for k — 3,4 the coefficients are given by the 
same expression (I19f) . The packing fraction at the I-N t 
bifurcation point, obtained as the solution of J23J, is 



In a second step we have calculated which gives 



1G 



ii 



x\x* 2 (y*)' 



■S3 



' (cosh Jfc 9) 
(cosh 



k fl\2 



-1, k = 1,2, 



11 (cosh 



2 ((^ + ct) 2 ) 

157T (w) 



16 (cosh 4 6»)(cosh^ 0) - 5 (cosh 3 0) 2 ' 



(44) 
(45) 
(46) 



(40) 



while the value corresponding to the isotropic-uniaxial 
nematic (N u ) bifurcation point, calculated directly from 



The region in the plane where H^ t < defines 

the demixing region, bounded by a dashed line in the 
inset of Fig. E^a) for k% = 1. As can be seen this region 
is closed as in the case of HDR for the FN transition. 
The main reason for this behavior is that in the one- 
component limit the FN t transition is of second order 



> and 



= 1 



2 {{L-af) 

37T (v) 



> 



(41) 



for all «i) which in 
do not 



(^) Nt 

turn means that the solutions of (x\ , K\ 
intersect the lines x* = {0, 1}. 

As the third and last step we have calculated the 
expressions for B^ and , which are too large to 
The equality of (H0J) and (gTJ defines a line in the plane be shown here, except for their one-component limits 
x\-K\ (fixing K2) where the FN t transition preempts the x* — > 1, 
I-N u transition. Taking the areas of both particles to be 
unity, as was done for HDR, and changing to variables 
9 V = ln(K„)/2 (v = 1,2), we obtain the solution of f?N t — 
?/n u as 



lim Bn u = 



1 5-2cotlr 6>i 
64 5-coth 2 (9i ' 



(47) 



15 coth 2 6»i-6 coth 4 6>i - 5 



2 cotlr 



1 -4sinh 2 6 2 
4 (sinh 2 6 1 - sinh 2 6» 2 ) ' 



(42) 



(48) 



with the constraint sinh 2 1 
v). 



< 1/4 < sinh 2 6V (and fi ^ 



This result should be taken with some care because 
the I-N Q (a = u, t) transitions can be of first order and 
the relative position of binodals can change the scenario 
predicted above. In order to elucidate the nature of these 
transitions, we have first calculated the coefficient B^ t for 
the tetratic phase which results in 



where, rather than H, the inverse isothermal compress- 
ibility zc^r was chosen as the adequate thermodynamic 
variable to elucidate the nature of the I-N„ transition of 
a one-component fluid |36| . An important difference to 
be mentioned is that, for HR, in the one-component 
limit can be negative for some n\ [see Eq. I|47f) ] , which in- 
dicates the presence of a tricritical point. In this limit the 



condition 



(xr T 



= is more stringent to determine 

the exact location of this point which gives, equating the 
numerator of the right-hand side of Eq. (|48|l to zero, the 
solution 



[22(11 + 5zi)r 2 - lQ Zl z 2 r + zx(U + 5z 2 )l 

4°zi 



(43) 



+ v / 7jl5 + 2/V5 



5.44 



where r = 
x* cosh 2 6> M /(cosh 2 6>) 



cosh 82 1 cosh 61 and 



Again the second-order poly- 
nomial with respect to r enclosed by the square brackets 
is always greater than zero. Thus B* > and the N t 
energy branch bifurcates from below from the I branch. 



[34j. On the other hand, under the assumption that 
both the I-N M and I-N t are of second order, the tran- 
sition to the N t phase preempts the FN M transition for 
Ki less than a value k* obtained from Eq. (|42|) by 
taking the one-component limit x\ = 1. The result is 
k* = (3 + \/5) /2 w 2.62. However, for this value of «*, 
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FIG. 3: Demixing regions predicted from H^ a = (a = u, t) 
for HR in the x\ — k,\ plane. The solutions of the equations 
H£ = (solid line), = (dashed line), and the function 
1421 1 (dotted lines) are shown for K2 = 1 (a), 5 (b), and 10 (c). 



the I-N„ transition is of first order, so that the value of 
k\ that determines when the N t phase begins to be sta- 
ble (for K\ < hi\) should be found as the intersection of 
the I-N t spinodal with the I binodal of the I-N„ coexis- 
tence. This value turns out to be [33| K i = 2.21. Thus, 
in the one-component limit, the range where the transi- 
tion from the I phase to an orientational ordered phase 
is of first order is 2.21 < ki < 5.44. For the mixture, the 
more stringent condition turns out to be H£ (x\, Ki) = 
whose solutions are shown for K2 = 1, 5, and 10 (see 
Fig. 0). The main differences between these figures and 
that obtained for HDR can be summarized as follows: (i) 
Some of the predicted demixing regions shown in Fig. |3| 
are open due to the first order nature of the I-N„ transi- 
tion in the one-component limit, as was discussed above. 
In contrast the solution to H^ t = generates a closed 
loop due to the second order nature of the I-Nj transi- 
tion [see the inset of Fig. Eta)]. Finally for k = 10 [Fig. 
0{c)] i/^ u = has two separate solutions, one of which 
is closed, the other one bounding an open region. In the 
Onsager limit we again obtain the closed loop shown in 
Fig. ^ as should be expected, since both models have 
the same asymptotic limit, (ii) Due to the presence of 
the N t phase in the mixture of HR, one of the curves 
which bounds the demixing regions is given by Eq. I|42|) 
(the dotted lines in Fig. (iii) The demixing regions 
for HR mixture are in general wider as compared to those 
of HDR. 

Let us discuss the behavior of the pressure, using Fig. 
|3fa) and, as an example, the case K2 = 1 and ni = 10 2 . 
This behavior is similar to that shown in Fig. [2] for 
the HDR model, namely the pressure is a monotonically 
decreasing function of composition. Also there are two 
points obtained from the intersection of the vertical line 
at K^ 1 = 10~ 2 and the boundaries of the shaded region 
(the upper one on the solid line and the other on the lower 
branch of the dashed line). These in principle would cor- 
respond to tricritical points as in the case of HDR. The 
point corresponding to a lower pressure (higher compo- 
sition) is a genuine tricritical point while the other one 
might be inside I-N M or I-N t demixing regions. 



C. The Zwanzig model 

This simple model, as applied to HR, allows the cal- 
culation of the phase diagram without any parametriza- 
tion because the orientational distribution function can 
be taken as h^) = [(1 + g„) /2] 5(4>) + [(l - g„) /2] 5{<j>- 
7r/2) \8{x) is the Dirac delta function], corresponding to 
a binary mixture of two species with perpendicular ori- 
entations. Thus, the excess part of the free energy is 
a second degree polynomial of the order parameters q v 
(with —1 < q v < 1), and the minimization of the total 
free energy requires to solve two transcendental equations 
to find their equilibrium values. Since no parametriza- 
tion is necessary, the location of all tricritical points, as 
obtained from the minimization, are exactly the same as 
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those obtained from the bifurcation analysis. This in fact 
is the reason why we have chosen to use this model. 

Once we calculate the phase diagram we can compare 
the results to those obtained using the bifurcation anal- 
ysis. To implement the latter we need the following ex- 
pressions: 



i« = 

,.( 2 ) = 



(<L„ - 2yx„ sinh 6» M sinh 9 V ) , 



a (3) =0 a 



(4) _ 



12' 



(49) 
(50) 



which allows, using the same procedure described above, 
to find the packing fraction at the bifurcation point, 



if = (1 + 2 (sinh 2 (9)) 1 , 



and the coefficient B* 



B* = 



1 (sinh 4 6») 
12 sinh 4 0, 



>0, 



(51) 



(52) 



where the constraint of all particle areas being equal 
to 1 was imposed and the same change of variables 
[Ofj, = ln(K M )/2] was used. Finally, to find the demixing 
behavior of this model, we need to make the expression 



U + y*) 4 



x\x* 2 (y*) 
x [l + (as 2 f + a 



l + a-|(l + a§) 



(as 2 - si) 2 



Sk 



S3 



(cosh 



2k , 



(cosh fe 6») 2 

' (sinh 4 9) 
(sinh 2 (9) 2 



1, fc=l,2 



1, a 



(cosh 2 6) 
(sinh 2 9) 



(53) 
(54) 
(55) 



to vanish. The solutions of this equation for different val- 
ues of K2 are shown in Fig. 21 As can be seen the general 
topology is similar to that found for HDR. Comparing 
both figures (^and^J we can draw as a conclusion that 
the discretization of orientations unfavors the FN demix- 
ing. 

The phase diagram of a Zwanzig binary mixture of 
species with aspect ratios k m equal to 9 (/x = 1) and 2 
(fx = 2) is shown in Fig. [S] in the pressure-composition 
plane. The tricritical point (where the I-Ni demixing 
transition changes from first to second order as the pres- 
sure is reduced) at x* = 0.1034 and (3P*v 1 = 5.2928 
coincides exactly with the point predicted from the bi- 
furcation analysis (see in Fig. 0] the upper intersection 
of the curve for k,i = 2 with a vertical line located at 
A = 2/9 w 0.222). 

An important feature of this phase diagram is that 
the binodals of the I-Ni or N1-N2 transitions tend to 
the one-component asymptotes x\ = and xx = 1 as 
the pressure is increased. This trend, characteristic of 
demixing scenarios (such as the usual I-I fluid demixing) , 
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FIG. 4: Demixing regions (shaded areas) of the Zwanzig 
model predicted from = 0. Different curves are the solu- 
tions of 1531 for K2 = 1, 2, 3, 5, and 00 (the Onsager limit) 
from outside to inside. 
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FIG. 5: Phase diagram of the Zwanzig mixture with m = 9 
and K2 — 2. 



should be compared with other possible topologies of 
phase diagrams (for example, those of three-dimensional 
mixtures of particles of similar lengths), where the first- 
order FN transition has two binodals that meet at the 
one-component limits of the phase diagram at finite pres- 
sure. This criterion is the one we have used to consider 
the transitions found here as FN demixing transitions 
instead of the standard FN orientational transition. 
As can be seen from Fig. \5\ at higher pressure from 



10 



the tricritical point there is a triple intersection between 
two first order transition lines (the I binodal of the I-Ni 
transition and the N2 binodal of the N1-N2 transition) 
and the second order I-N2 transition line. This point is 
called in the literature a critical endpoint. An important 
remark to make is that this critical endpoint located at 
x* = 0.0144 and (3P*vx = 16.2775 does not exactly coin- 
cide with the predicted result from the bifurcation analy- 
sis. If we look at the lower intersection between the curve 
corresponding to k 2 = 2 and the vertical line located at 
A 0.222 in Fig. 0] we find the value x* which corre- 
sponds to PPvi — 13.8014, which is lower than the above 
mentioned critical endpoint pressure. The reason for this 
disagreement can be elucidated if we plot the Gibbs free 
energy per particle of the mixture as a function of x\ for 
three different values of the pressure, f3Pv\ — 13.8014, 
which is the pressure corresponding to the upper tricrit- 
ical point as obtained from the bifurcation analysis [the 
open circle in the inset of Fig. |Sf a) corresponds to the bi- 
furcation point]; PPvi — 16.2775, which corresponds to 
the triple intersection [Fig. EJb)] where the bifurcation 
and the coexistence points coalesce; and (3Pv\ — 18.7536, 
where the two coexistence points already lie on the ne- 
matic branch. 

The main conclusion we can draw from this is that the 
bifurcation analysis predicts exactly the tricritical point 
(which corresponds to the lower pressure or higher x\ in 
Fig. [3} and gives an approximate value for the critical 
endpoint. 

In the light of these results we can interpret Fig. 02a) 
and try to predict, as an example, the demixing behav- 
ior of a mixture of freely rotating HR with K2 = 1 and 
Ki = 10 2 . If we increase the pressure j3Pv\(x\) following 
the I-N„ spinodal curve we find that at a composition 
defined as the intersection with the solid line (the tri- 
critical point) the mixture begins to demix into I and 
N M phases. The demixing gap becomes wider with pres- 
sure, and at high pressures the lower intersection with the 
dashed line defines the approximate location of a critical 
endpoint where the first order I-N u and N u -N t transitions 
coalesce with a second order I-N t transition. Increasing 
further the pressure we should have a N u -N t coexistence 
and, ultimately, a coexistence between two uniaxial ne- 
matic phases. This latter coexistence could be predicted 
if we take into account the phase diagram for the one- 
component system, where the N t phase is sandwiched 
between the I and N„ phases [33|. 

V. PHASE DIAGRAMS OF HDR MIXTURES 

In this section we want to explicity discuss the demix- 
ing scenarios that may occur in binary mixtures of hard 
particles. We will focus our attention on HDR mixtures. 
For this purpose we have parameterized the orientational 
distribution functions of each species as 

M0) = —n — \ ex P (°> cos 2< /0 i ( 56 ) 
t^o ("/J 
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FIG. 6: Gibbs free energy per particle in reduced units /3g 
as a function of x\ for the Zwanzig mixture with species 
having ki = 9 and K2 = 2. The values of the pres- 
sure are /3Pvi = 13.8014 (a), 0Pvi = 16.2775 (b), and 
j3Pv\ = 18.7536 (c). The solid and dotted lines correspond 
to the I and N branches, respectively. The filled circles rep- 
resent the coexistence points. The insets are enlargements of 
the neighborhood of the bifurcation points, which are shown 
by open circles. Panel (a) shows the loss of convexity of j3g 
in the N brach at the bifurcation point. 
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FIG. 7: (a) Phase diagram of the mixture of hard disks k,2 = 
1 and HDR with ki = 10. (b) A zoom taken around the 
tricritical point. 
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where Iq {x) is the zeroth-order modified Bessel function. 
This parametrization fulfills the normalization constraint 
fo d4>h^4>) = 1. The free energy per volume (see Sec. 
II for its expression) was minimized with respect to 
(ji = 1,2) and the double-tangent construction on the 
thermodynamic potential (3g(x\) (the Gibbs energy per 
particle in reduced units as a function of the composition 
of the mixture for a fixed pressure) was used to calculate 
the coexistence values of the composition of the mixture 
and the packing fraction. 

The results for a mixture of hard disks (k2 = 1) and 
HDR with Ki = 10 are shown in Figs. |7£a) and[7Jb). Fol- 
low pressures there is a continuous I-N transition at a 
pressure which agrees with that calculated from the spin- 
odal packing fraction curve (|34[1 . At some composition 
(x\ = 0.393) we find a tricritical point from which the 



FIG. 9: Phase diagram for a HDR mixture with ki = 4 and 

K2 = 2. 



mixture begins to demix into the I phase (rich in discs) 
and the N phase (rich in rods). The location of the tricrit- 
ical point should be compared with that resulting from 
the bifurcation analysis (x* = 0.434) which gives the ex- 
act result. The difference is due to the parametrization 
used. As we have already pointed out the value of x\ de- 
pends on the coefficient B* which in turn depends on the 
d$ (k — 1, . . . ,4). These are the expansion coefficients 
of the free-energy difference around the bifurcation point 
up to fourth order. Although the used parametrization 
captures the right second order terms, further terms are 
only approximate due to the restriction on the minimiza- 
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tion variables to be only one per species. 

An interesting result is shown in Fig. Efb), a zoom of 
Fig. Ufa) around the tricritical point. As we can see this 
part of the phase diagram exhibits a N1-N2 coexistence 
which ends in a critical point, along with the presence of 
a I-N1-N2 triple point. As these features are similar to 
those observed in the Zwanzig mixture (not shown here) , 
where no parametrization was used, we are confident that 
this scenario is qualitatively correct. 

Two other phase diagrams are shown in Figs. [S]and 
for different mixtures. The topology of the first one 
is similar to that found in the Zwanzig mixture, where 
the I-Ni demixing is followed by N1-N2 coexistence and 
the phase diagram includes one tricritical and one critical 
endpoint. 

The third phase diagram (see Fig. [3J includes an I-N 
transition which is always of second order and, at high 
pressures, a N1-N2 demixing ending in a critical point. 
This kind of phase diagram is typical of mixtures of par- 
ticles with similar aspect ratios. Since the pressures at 
which the demixing transition occurs are rather high, 
we expect this transition to be metastable with respect 
to phases with partial or full spatial order. Comparing 
the phase diagram of Fig. |5| with the exact bifurcation 
analysis results (see Fig. H), we can conclude that, even 
though I-N demixing is already allowed from A < 0.595, 
the parametrization used changes the demixing behavior 
of the mixture. An exact free-energy minimization for 
the same value of K2 and for A > 0.595 should qualita- 
tively give a phase diagram similar to that shown in Fig. 

El 

VI. CONCLUSIONS 

While I-I demixing is forbidden in two-dimensional 
hard-body additive mixtures of anisotropic particles ac- 



cording to SPT 25] , we have shown in the present work 
that the inclusion of ordered phases with orientational 
symmetry breaking changes completely the demixing sce- 
nario. The phase diagrams characteristic of these mix- 
tures can exhibit I-N and N-N demixing. We have used 
a bifurcation analysis to demonstrate rigorously that I-N 
demixing occurs. On the other hand, explicit calcula- 
tions of the phase diagrams of HDR, using an accurate 
parametrization, have been performed in order to show 
the occurrence of N-N demixing. The simple structure 
of SPT allowed us to obtain analytically the stability cri- 
terion for the mixture as a function of packing fraction, 
composition and shape of the constituent particles. Thus, 
using this procedure, we can predict the I-N demixing 
scenarios for mixtures of HDR and HR. 

To show the relative stability of these demixed phases 
with respect to nonuniform phases (e.g., solid phases), we 
would need to carry out a full minimization, with respect 
to the density profile /fy(r, <f>) (which also depends on 
spatial variables), of a density functional constructed in 
such a way that it recovers the SPT in the uniform limit. 
This work is a task in progress. 

Finally, we expect that all the demixing scenarios pre- 
dicted here be confirmed by computer simulations of 
hard-body mixtures in two dimensions. These simula- 
tions are still lacking. 
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